clear
set more off
local n = 2
local k = 0.8
// author: Lucas

cap log close
log using "$logdir/analysis/progressivity.log", replace

local incvar soi_agi

use $statadir/soi_cdw_taxsim_2019 if nonfiler == 0, clear

// calculate SOI tax liability

*****************************************
*										*
*				Tax errors				*
*										*
*****************************************
gen double delta = taxsim1_fed_incm_tax - soi_fed_incm_tax
qui drop if missing(delta) // bad tax unit matches
qui drop if soi_agi <= 0

fasterxtile bin = soi_agi  [aw=soi_wgt], nq(100)

qui replace bin = ceil(bin/5)*5 if bin <= 95


/*
// disclosure
gen abs_delta = abs(delta)
bys bin: gegen tot_delta = sum(abs_delta) [iw=soi_wgt]
gen sortvar = -abs_delta
sort bin sortvar
by bin: gen cumul = sum(abs_delta*soi_wgt)/tot_delta
by bin: gen bad = _n == `n' & cumul >= `k'
qui sum bad
assert r(max) == 0
*/

collapse m_agi = soi_agi (sum) sim = taxsim1_fed_incm_tax true = soi_fed_incm_tax delta soi_agi [iw=soi_wgt], by(bin) fast
foreach var in sim true delta soi_agi {
	qui replace `var' = `var'*10^-9
}

gen share = delta/abs(soi_agi)
list bin share m_agi
keep if bin > 1
qui replace m_agi = m_agi/1000
niceloglabels m_agi, style(13) local(labname)

twoway scatter share m_agi, xtitle("Income (000s)") ytitle("Pre-population error as share of income") ///
	scheme(s1color) xscale(log) xlabel(`labname') yline(0, lpattern("--"))
graph export $outdir/figures/progressivity.pdf, replace
graph export $outdir/figures/progressivity.eps, replace

cap log close